﻿using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Linq;
using System.Numerics;
using System.Text;
using System.Threading.Tasks;
using MathNet.Numerics.IntegralTransforms;
using phminarm.algorithm.Base;

namespace phminarm.algorithm.DataProcess
{
    public class FFTWraper
    {
        public static float recover { get; private set; } = 2;

        public static (float[], float[]) GetFFTAmp(float[] timewave, WindowType windowType = WindowType.Hanning)
        {
            Stopwatch sw = new Stopwatch();
            sw.Start();

            int N = timewave.Length;
            //N不是2^n的情况处理

            //数据格式转换
            Complex[] data = new Complex[N];

            //加窗
            float[] inputW = Windowed(timewave, windowType);
            for (int i = 0; i < N; i++)
            {
                data[i] = inputW[i];
            }

            //快速傅里叶变换
            //using (var pinIn = new PinnedArray<Complex>(input))
            //using (var pinOut = new PinnedArray<Complex>(output))
            //{
            //    //DFT.Wisdom.Import("wisdom.txt");
            //    DFT.FFT(pinIn, pinOut);
            //}
            Fourier.Forward(data, FourierOptions.Matlab);

            //将双边谱转化为单边谱,汉宁恢复系数为2
            switch (windowType)
            {
                case WindowType.Hanning:
                    recover = 2;
                    break;

            }
            float[] spec = data.Select(x => Convert.ToSingle(x.Magnitude / (N / 2.0) * recover)).Take(N / 2).ToArray();
            spec[0] = spec[0] / 2;
            float[] pha = data.Select(x => Convert.ToSingle(x.Phase)).Take(N / 2).ToArray();

            sw.Stop();
            TimeSpan ts2 = sw.Elapsed;
            Console.WriteLine("快速傅里叶变换总共花费{0}ms.", ts2.TotalMilliseconds);
            return (spec, pha);
        }

        public static (float[], float[]) GetFFTPower(float[] timewave, WindowType windowType = WindowType.Hanning)
        {
            int N = timewave.Length;
            //N不是2^n的情况处理

            //数据格式转换
            Complex[] data = new Complex[N];

            //加窗
            float[] inputW = Windowed(timewave, windowType);
            for (int i = 0; i < N; i++)
            {
                data[i] = inputW[i];
            }

            Fourier.Forward(data, FourierOptions.Matlab);

            //将双边谱转化为单边谱,汉宁恢复能量系数为1.633
            float recover = 0;
            switch (windowType)
            {
                case WindowType.Hanning:
                    recover = 1.633f;
                    break;

            }
            float[] spec = data.Select(x => Convert.ToSingle(x.Magnitude / (N / 2.0) * recover)).Take(N / 2).ToArray();
            spec[0] = spec[0] / 2;
            float[] pha = data.Select(x => Convert.ToSingle(x.Phase)).Take(N / 2).ToArray();

            return (spec, pha);
        }

        public static void Test()
        {


        }

        public static float[] Windowed(float[] input, WindowType windowType = WindowType.Hanning)
        {
            int N = input.Length;
            float[] output = new float[N];
            switch (windowType)
            {
                case WindowType.Hanning:
                    for (int i = 0; i < N; i++)
                    {
                        output[i] = Convert.ToSingle(input[i] * (0.5 + 0.5 * Math.Cos(2 * Math.PI * i / (N - 1))));
                    }
                    break;
                case WindowType.Uniform:

                    break;
            }
            return output;
        }
    }
}
